Renormalization of earthquake aftershocks 

Anne Sornette^ and Didier Sornette^'^ 
^ Lahoratoire de Physique de la Matiere Condensee, CNRS UMR 6622 
Universite de Nice-Sophia Antipolis, Pare Valrose, 06108 Niee, France 
^ Department of Earth and Space Sciences and Institute of Geophysics and Planetary Physics 
University of California, Los Angeles, California 90095-1567 

Abstract 



Assume that each earthquake can produce a series of aftershock indepen- 
dently of its size according to its "local" Omori's law with exponent l-\-9. Each 
aftershock can itself trigger other aftershocks and so on. The global observable 
Omori's law is found to have two distinct power law regimes, the first one with 
exponent p- = 1 — 6 for time t < t* ~ k~^^^ , where < 1 — k < 1 measures the 
fraction of triggered earthquakes per triggering earthquake, and the second one 
with exponent p+ = 1 + for larger times. The existence of these two regimes 
rationalizes the observation of Kisslinger and Jones [1991] that the exponent 
p seems positively correlated to the surface heat flow: a higher heat flow is a 
signature of a higher crustal temperature, which leads to larger strain relax- 
ation by creep, corresponding to fewer events triggered per earthquake, i.e. to 
a larger k, and thus to a smaller t*, leading to an effective measured exponent 
more heavily weighted toward p+ > 1. 
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I. INTRODUCTION 



Aftershocks constitute a very significant percentage of even the large shallow earthquakes 
[Hough and Jones, 1997] and all the more so for the small earthquakes that dominate the 
catalogs. An aftershock is usually defined as any earthquake that occurs within a distance 
equal to the length of the fault that ruptured during the mainshock and during the span of 
time that seismicity rate in that region remains above its pre- mainshock background level. 

This definition can in fact be imprecise for a variety of reasons : 

1. the pre-mainshock rate can be ill-defined due to inherent intermittent behavior [Kagan 
and Jackson, 1991; Kagan, 1994] or poorly constrained by sparse data. The intermit- 
tency or clustering in time is sometimes quantified by fractal tools [Dattatrayam and 
Kamble, 1994; Henderson, 1994]. 

2. Some earthquakes are sometimes triggered hundreds of kilometers away, at significantly 
greater distance than the mainshock rupture length [Hill et al., 1993; Steeples and 
Steeples, 1996]. 

3. Patterns of spatial and temporal self-organization between events at large distance and 
time scales have been observed [Knopoff et al., 1997; Bowman et al., 1998] that suggest 
more complex interactions than a simple time delay mechanism between mainshock 
and aftershocks, and has been ascribed to a kind of self-organized critical behavior 
[Sornette and Sornette, 1989; Bak and Tang, 1989; Grasso and Sornette, 1998]; 

4. In some models of self-organized spatio-temporal earthquake behavior [Huang et al., 
1998; Hainz et al, 1998], the same physical mechanism is found to produce both the 
Gutenberg- Richter frequency-magnitude distribution and the foreshock and aftershock 
Omori's laws. 

The identification of an aftershock becomes feasible, by definition, only after an event has 
been recognized as a main shock. This main shock can itself be an aftershock of a larger event 



that preceded it. The time-delayed interactions between mainshock and aftershocks, that 
are usually invoked [Scholz, 1990; Das and Scholz, 1981; Shaw, 1993], have no reason not to 
be at work between two events when the second one is the large. Kagan and Knopoff [1978] 
have even asserted that both foreshocks and aftershocks are a manifestation of essentially 
the same process. 

These considerations lead us to formulate the following hypothesis and explore its self- 
consistency : 

Hypothesis : extending the usual definition of an aftershock to allow for the fact that the 
delayed triggered events can be of arbitrary size and not just smaller, can we construct a 
coherent model of aftershocks in which each earthquake is the generator of a sequence of de- 
layed triggered events of arbitrary size, in other words, each earthquake is the "generalized" 
aftershock of a some preceding event? 

To address this question, we explore the following model. Each event triggers its own 
"local" Omori's rate of aftershocks. Each of these aftershocks itself triggers its own "local" 
Omori's law of aftershock, and so on in an infinite branching process. The relevant question 
for observations is what is the resulting global Omori's law? Indeed, the individual "local" 
Omori's rate triggered by each event is not observable, only the superposition of all events 
at all levels of the cascade is usually quantified. This problem belongs to the general class 
in which one specifies the laws at the "microscopic" scale, here the local Omori's law, and 
strives to derive the resulting "macroscopic" laws. In this paper, we do not address the 
mechanism(s) at the origin of the "local" Omori's law but just take it for granted. We 
rather address the collective behavior, i.e. how the "local" law is "renormalized" [Wilson, 
1979] into a macroscopic observable Omori's law. 
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II. SELF-CONSISTENT OMORI'S LAW 



To formulate the question in formal terms, we follow Kagan and Knopoff (1981) and 
assume that the probability that one event occurring at time zero gives birth to another in 
the time interval between t and t + dt is 

^(t)dt=il-K)9 4-^H{t-to), (1) 

where H{t) is the Heaviside function: H{t — to) = for t < to and 1 otherwise, which means 
that the triggering of seismic activity starts a short delay to after the shock. This minimum 
waiting time is necessary for the total number of events to be finite, otherwise 0(t) would 
not be normalized due to the divergence at t — > 0. Physically, it represents the simplest way 
to account for a progressive build up of activity after the main shock. 

This expression is exactly Omori's law for the rate of aftershocks following a main shock, 
albeit with the modification that we do not specify that aftershocks have to be smaller. The 
exponent 1 + 6 of the "local" Omori's law has not reason a priori to be the same as the one 
measured macroscopically which is usually found between 0.8 and 1.2 with an often quoted 
value 1. This is in fact the question we address: assuming the form (P for the "local" 
Omori's law, is the global Omori's law still a power law and, if yes, how does its exponent 
depend on 6^? 

The integral of 0(t) over time is the average number of earthquakes created by each event 
and is equal to 

r+oo 

/ (f){t) dt = l-K. (2) 
Jo 

If K < 0, more than one earthquake is triggered per earthquake. This regime corresponds to 
the super-critical regime of branching processes [Harris, 1963] in which the total number of 
events grows exponentially with time. If k > 0, there is less than one earthquake triggered 
per earthquake. This is the sub-critical regime in which the number of events following 
the first main shock decays eventually to zero. The critical case k = is at the borderline 



between the two regimes. In this case, there is exactly one earthquake on average triggered 
per earthquake and the process is exactly at the critical point between death on the long 
run and exponential proliferation. In the sequel, we take k > as the physically relevant 
regime where the aftershock activity following a main shock eventually decays at long times. 

We analyze the case where there is an origin of time t = at which we start recording 
the rate of earthquakes, assuming that the largest earthquake of all has just occurred at 
t = and somehow reset the clock. In the following calculation, we will forget about the 
effect of events preceding the one at t = and count aftershocks that are created only by 
this main shock. 

Let us call N{t)dt the number of earthquakes between t and t + dt. N{t) is the solution of 
a self-consistency equation that formalizes mathematically the following process : the first 
earthquake may trigger aftershocks; these aftershocks may trigger their own aftershocks, 
and so on. The rate of seismicity at a given time t is the result of this cascade process. The 
self-consistency equation that sums up this cascade reads 



Its meaning is the following. The rate N{t) at time t is the sum over all induced rates from 
all earthquakes that occurred at all previous times with an induced rate per earthquake 
occurring at time t — t equal to 0(t — r). The essential point is that the same function 
appears on the l.h.s. and in the integrant in the r.h.s. 

The lower bound in the integral in (^) expresses the fact that we look at the time 
development of the seismic activity after an origin of time, that can be taken for instance as 
the triggering time at which a big event occurred. Such an origin of time is important as we 
are not describing a stationary process, but rather a rate which on average relaxes (slowly) 
to zero. 

The problem is then to determine the functional form of N{t), assuming that (p is given 
by The integral equation is a Wiener- Hopf integral equation [Feller, 1971]. It is 




(3) 
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well-known [Feller, 1971; Morse and Feshback, 1953] that, if 0(r) decays no slower than an 
exponential, then N{t) has an exponential tail N{t) ~ exp[— /it] for large t with /i solution 
of / 0(x) exp[yux] dx = 1. In the present case, 0(r) decays much slower than an exponential 
and a different analysis is called for that we now present. 

III. SOLUTION OF THE WIENER-HOPF INTEGRAL EQUATION FOR POWER 

LAW KERNEL 

Inserting (|I|) in (^) gives 

N{t) = (1 - k) N{r) etUt- t)-'^^+'^ dr. (4) 

The upper bound t — to comes from the Heaviside function and reflects the rule that the 
rate N{t) at t is influenced only by events that occurred at least a time to before. 
By writing (H) for time t + "we get 

N{t + to) = (1 - /t) f N{t) Kt,{t - r)dT , (5) 
Jo 

where we define the Kernel 

Kt,{x)^d4ix + to)-'^'+'K (6) 

We recognize that the r.h.s. of (|^) is exactly the convolution of N{t) with Kt^Xt)- Since 
we are dealing with a non- stationary process (i.e. there is an origin of time) and we have 
a convolution operator, the natural tool is the Laplace transform /(/3) = J^"^ f{t)e~^^dt. 
Applying the Laplace transform to (|^) yields 

/ N{t + to)e-^'dt = (1 - K)N{f3)Kt,{(3) , (7) 
Jo 

where the r.h.s. has used the convolution theorem. To deal with the l.h.s., we transform it 

as 



r+oo / r+oo rto \ 

N{t + to)e~'^'dt = e^^° / N{t)e-^^dt = e^*« / N{t)e-^^dt - / N{t)e-^^dt . 

Jto \Jo Jo J 

(8) 

Recall that the first great event is assumed to occur at t = 0. During the time interval from 
to to, no other event occurs. Thus, /q" N{t)e'^Ht = toN{0) with A^(0) = 1. We also make 
the approximation e'^*° ~ 1, based on the condition that we look at earthquake rates at 
long times to << t which translates into the condition /? << tg ^ for the Laplace conjugate 
variable p. We thus get 

N{t + to)e-^'dt ^ N{p) - toN{0) . (9) 
Reporting into and solving for N{i3) gives the expression of the Laplace transform of 

N{ty. 

m ^ '-m . (10) 

The solution for N{t) is derived by taking the inverse Laplace transform of (p^). 

In this goal, we shall now make use of a well-known result for the expression of the 
Laplace transform of a power law. The important point to notice is that K(t) ~ for 
large t. Then, its Laplace transform takes the form (obtained by integrating by part / times, 
where / the integer part of 6* i.e. I < 9 < I + 1) 

I ^-i^-^(^-i)(^^-2)...(^-/)y^ 

+ — ^-'^yf^' r dxe-^x'-r (11) 



-i)[e-2)...{e-i) J(3 

This last integral is equal to 

(3^ dxe-'^x^-^ = T {I + 1-6)1(3^ + (3^+^ -f* {I + 1-6,(3)], (12) 
where F is the Gamma function (r(n + 1) = n!) and 



+00 on 

is the incomplete Gamma function. We see that K{P) presents a regular Taylor expansion 
in powers of (3 up to the order /, followed by a term of the form f3^. We can thus write 

k{f3) = l+ri/3 + + ri(3' + ref3'> + (14) 



with ri = — (t), r2 = ... are the moments of K. For small (3, we rewrite K{(3) under the 
form 



K{p) = exp 



(15) 



Lk=l 

where the coefficient can be simply expressed in terms of the r^'s. We recognize that 
this transformation is similar to that from the moments to the cumulants. The expression 
([T5|) generalizes the canonical form of the characteristic function of the stable Levy laws 
[Gnedenko and Kolmogorov, 1954] for arbitrary values of 6, and not solely for ^ < 2 for 
which they are defined. The canonical form is recovered for 6 < 2 for which the coefficient 
(^2 is not defined (the variance does not exist) and the only analytical term is {t)(3 (for 
/i > 1). 

Here, we are interested in the case < 6' < 1 for which 

KtoiP) = I - diPtoY + h.o.t. (16) 
where h.o.t. stands for higher order terms. Reporting in (p!OD this gives 

m = = . (17) 



Two cases must be distinguished. 

- (iStoY « K, corresponds to t >> t* = \. In this case, we can expand , 



k9 

which leads to 



iV(/3) = ^[l-rfi^(/3to)% (18) 



We recognize the Laplace transform ([15|) of a power law of exponent 6', i.e. 

N{t) ~ r(^+^) for t»t* . (19) 



For t < t* , d{l — K)(pto) > K and ([T7| ) can be written with a good approximation as 



m = ''""[Z - iPtor. (20) 

Denoting r{z) = J^l''^ dt e"* we see that /o^°° cit e"^* t^-^ = r(z)/3-^ Comparing with 
), we thus get 

N{t) ~ r^^-^^ for t < r . (21) 



IV. DISCUSSION 

Treating all aftershocks on the same footing by assuming that each may trigger new 
aftershocks with the same "local" rate gives a global Omori's law N{t) described by two 
power laws with exponents = 1 — 6 for t < t* and = 1 + 6 for t > t* , where the 
characteristic time t* = t^/K}/^ . 

• Our results shows that the simple and parsimonious model in which all earthquakes are 
put on the same footing, i.e. are all susceptible to trigger its train of aftershocks is fully 
consistent with empirical observations. It suggests that the usual taxonomy identifying 
aftershocks as special events may not resist a more physically-based analysis. 

• The value 0^0, corresponding to a local Omori's exponent equal to one, is the only 
value that gives a completely self-similar rate, i.e. the same power law for all times, and in 
addition the same decay rate at the "local" and global level. In other words, the 1/t law is 
the only one that solves the self-consistent condition that the same law describes all levels 
of aftershock triggering processes. 



• In practice, even a large deviation from ^ = 0, say, 6 — 1/2 may produce a reasonable 

approximation to the usual 1/t Omori's power law, due to the existence of the long cross-over 
between the two power laws whose exponents possess the value 1 as their exact barycenter. 
In this picture, the observed Omori's law is not exactly 1/t but a mixture of two power 
laws slowly crossing over from one to the other. This can also provide a mechanism for the 
variations of measured Omori's exponent [Kisslinger and Jones, 1991]. 

• Let us take 9 — 0.2 and k — 10~^. This yields t* — 10^ for the characteristic time scale. 
If to = 1 s, then t* ^ 1 day. Suppose now that n = 0.05. This yields t* ^ 3.2 10^ to ^ 30 
days. Thus, a small variation of the fraction 1 — k of earthquakes triggered per earthquake 
(from 0.9 to 0.95) leads to a dramatic variation (day to month) of the characteristic time 
scale t*. This offers a simple and physically appealing explanation for the observations of 
Kisslinger and Jones [1991] that Omori's exponent p seems to be positively correlated to 
the surface heat flow : since a higher heat flow is probably a signature of a higher crustal 
temperature, this leads to larger strain relaxation by creep. As a consequence, less strain 
is accomodated by earthquakes and there are thus fewer events triggered per earthquake, 
meaning that k is larger and t* is smaller. Now, the effective exponent p measured over a 
time interval t is approximately an average of p_ and p+ weighted respectively by t* — to 
and t — t*. If t* decreases (resp. increases), the effective exponent p increases toward p+ 
(resp. decreases toward p_). Our theory thus views Omori's exponent as a direct gauge of 
the fraction of strain accomodated by earthquakes. This prediction can be tested by future 
investigations. 

• Kagan and Knopoff [1981] have used a similar model to the one studied in the present 
paper in order to generate synthetic catalogs. Their model differs in the definition of the 
events, in the introduction of an additional exponential relaxation and in the introduction 
of a threshold for detection. These complications do not allow for an analytic treatment as 
performed here but their numerical simulations are compatible with our results, in particular 



with respect to the appearance of a long time scale t* — Iq/k}/^ . 



We acknowledge stimulating discussions with Y. Kagan and L. Knopoff. 
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